home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / lib / mathlib / libblas / src_original / dznrm2.f < prev    next >
Encoding:
Text File  |  1994-08-02  |  4.9 KB  |  169 lines
  1.       DOUBLE PRECISION FUNCTION DZNRM2( N, ZX, INCX )
  2. *
  3. *     unitary norm of the complex n-vector stored in zx() with storage
  4. *     increment incx .
  5. *     if    n .le. 0 return with result = 0.
  6. *     if n .ge. 1 then incx must be .ge. 1
  7. *
  8. *           c.l.lawson , 1978 jan 08
  9. *
  10. *     four phase method     using two built-in constants that are
  11. *     hopefully applicable to all machines.
  12. *         cutlo = maximum of  sqrt(u/eps)  over all known machines.
  13. *         cuthi = minimum of  sqrt(v)      over all known machines.
  14. *     where
  15. *         eps = smallest no. such that eps + 1. .gt. 1.
  16. *         u   = smallest positive no.   (underflow limit)
  17. *         v   = largest  no.            (overflow  limit)
  18. *
  19. *     brief outline of algorithm..
  20. *
  21. *     phase 1    scans zero components.
  22. *     move to phase 2 when a component is nonzero and .le. cutlo
  23. *     move to phase 3 when a component is .gt. cutlo
  24. *     move to phase 4 when a component is .ge. cuthi/m
  25. *     where m = n for x() real and m = 2*n for complex.
  26. *
  27. *     values for cutlo and cuthi..
  28. *     from the environmental parameters listed in the imsl converter
  29. *     document the limiting values are as follows..
  30. *     cutlo, s.p.   u/eps = 2**(-102) for  honeywell.  close seconds are
  31. *                   univac and dec at 2**(-103)
  32. *                   thus cutlo = 2**(-51) = 4.44089e-16
  33. *     cuthi, s.p.   v = 2**127 for univac, honeywell, and dec.
  34. *                   thus cuthi = 2**(63.5) = 1.30438e19
  35. *     cutlo, d.p.   u/eps = 2**(-67) for honeywell and dec.
  36. *                   thus cutlo = 2**(-33.5) = 8.23181d-11
  37. *     cuthi, d.p.   same as s.p.  cuthi = 1.30438d19
  38. *     data cutlo, cuthi / 8.232d-11,  1.304d19 /
  39. *     data cutlo, cuthi / 4.441e-16,  1.304e19 /
  40. *
  41. *     .. Scalar Arguments ..
  42.       INTEGER                           INCX, N
  43. *     ..
  44. *     .. Array Arguments ..
  45.       COMPLEX*16                        ZX( 1 )
  46. *     ..
  47. *     .. Local Scalars ..
  48.       LOGICAL                           IMAG, SCALE
  49.       INTEGER                           I, IX, NEXT, NN
  50.       DOUBLE PRECISION                  ABSX, CUTHI, CUTLO, HITEST, ONE,
  51.      $                                  SUM, XMAX, ZERO
  52.       COMPLEX*16                        ZDUMI, ZDUMR
  53. *     ..
  54. *     .. Intrinsic Functions ..
  55.       INTRINSIC                         DABS, DSQRT, FLOAT
  56. *     ..
  57. *     .. Statement Functions ..
  58.       DOUBLE PRECISION                  DIMAG, DREAL
  59. *     ..
  60. *     .. Statement Function definitions ..
  61.       DREAL( ZDUMR ) = ZDUMR
  62.       DIMAG( ZDUMI ) = ( 0.0D0, -1.0D0 )*ZDUMI
  63. *     ..
  64. *     .. Data statements ..
  65.       DATA                              ZERO, ONE / 0.0D0, 1.0D0 /
  66.       DATA                              CUTLO, CUTHI / 8.232D-11,
  67.      $                                  1.304D19 /
  68. *     ..
  69. *     .. Executable Statements ..
  70.       IF( N.GT.0 )
  71.      $   GO TO 10
  72.       DZNRM2 = ZERO
  73.       GO TO 140
  74. *
  75.    10 ASSIGN 20 TO NEXT
  76.       SUM = ZERO
  77.       IX = 1
  78.       IF( INCX.LT.0 )
  79.      $   IX = 1 - ( N-1 )*INCX
  80.       NN = IX + ( N-1 )*INCX
  81. *
  82. *        begin main loop
  83. *
  84.       DO 130 I = IX, NN, INCX
  85.          ABSX = DABS( DREAL( ZX( I ) ) )
  86.          IMAG = .FALSE.
  87.          GO TO NEXT( 20, 30, 60, 110, 70 )
  88.    20    IF( ABSX.GT.CUTLO )
  89.      $      GO TO 100
  90.          ASSIGN 30 TO NEXT
  91.          SCALE = .FALSE.
  92. *
  93. *           phase 1.  sum is zero
  94. *
  95.    30    IF( ABSX.EQ.ZERO )
  96.      $      GO TO 120
  97.          IF( ABSX.GT.CUTLO )
  98.      $      GO TO 100
  99. *
  100. *           prepare for phase 2.
  101. *
  102.          ASSIGN 60 TO NEXT
  103.          GO TO 50
  104. *
  105. *           prepare for phase 4.
  106. *
  107.    40    ASSIGN 70 TO NEXT
  108.          SUM = ( SUM / ABSX ) / ABSX
  109.    50    SCALE = .TRUE.
  110.          XMAX = ABSX
  111.          GO TO 80
  112. *
  113. *           phase 2.  sum is small.
  114. *                     scale to avoid destructive underflow.
  115. *
  116.    60    IF( ABSX.GT.CUTLO )
  117.      $      GO TO 90
  118. *
  119. *           common code for phases 2 and 4.
  120. *           in phase 4 sum is large.  scale to avoid overflow.
  121. *
  122.    70    IF( ABSX.LE.XMAX )
  123.      $      GO TO 80
  124.          SUM = ONE + SUM*( XMAX / ABSX )**2
  125.          XMAX = ABSX
  126.          GO TO 120
  127. *
  128.    80    SUM = SUM + ( ABSX / XMAX )**2
  129.          GO TO 120
  130. *
  131. *           prepare for phase 3.
  132. *
  133.    90    SUM = ( SUM*XMAX )*XMAX
  134. *
  135.   100    ASSIGN 110 TO NEXT
  136.          SCALE = .FALSE.
  137. *
  138. *           for real or d.p. set hitest = cuthi/n
  139. *           for complex      set hitest = cuthi/(2*n)
  140. *
  141.          HITEST = CUTHI / FLOAT( N )
  142. *
  143. *           phase 3.  sum is mid-range.  no scaling.
  144. *
  145.   110    IF( ABSX.GE.HITEST )
  146.      $      GO TO 40
  147.          SUM = SUM + ABSX**2
  148.   120    CONTINUE
  149. *
  150. *           control selection of real and imaginary parts.
  151. *
  152.          IF( IMAG )
  153.      $      GO TO 130
  154.          ABSX = DABS( DIMAG( ZX( I ) ) )
  155.          IMAG = .TRUE.
  156.          GO TO NEXT( 30, 60, 110, 70 )
  157. *
  158.   130 CONTINUE
  159. *
  160. *           end of main loop.
  161. *           compute square root and adjust for scaling.
  162. *
  163.       DZNRM2 = DSQRT( SUM )
  164.       IF( SCALE )
  165.      $   DZNRM2 = DZNRM2*XMAX
  166.   140 CONTINUE
  167.       RETURN
  168.       END
  169.